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Abstract 

We study the early time behavior of the Baxter- Wu model, an Ising model 
with three-spin interactions on a triangular lattice. Our estimates for the 
dynamic exponent z are compatible with results recently obtained for two 
models which belong to the same universality class of the Baxter- Wu model: 
the two-dimensional four-state Potts model and the Ising model with three- 
spin interactions in one direction. However, our estimates for the dynamic 
exponent 9 of the Baxter- Wu model are completely different from the values 
obtained for those models. This discrepancy could be related to the absence 
of a marginal operator in the Baxter- Wu model. 
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I. INTRODUCTION 



Since the work of Janssen et al. [1] and Huse [2] pointing out that universality and scaling 
are present at the early stage of the time evolution of dynamic systems, several statistical 
models have been investigated using Monte Carlo simulations in the short-time regime [3]. 
This kind of investigation avoids the critical slowing down effects and provides an alternate 
method for calculating the dynamic exponent z, the static exponents /3 and v as well as 
the new dynamic exponent 9 that governs the anomalous behavior of the nonequilibrium 
magnetization when the system is quenched to the critical temperature T c . Using this new 
approach kinetic Ising and Potts models [4-8] were exaustively studied in order to confirm 
available results and to enlarge our knowledge of critical phenomena, including the critical 
relaxation of the early time dynamics. 

In this paper, we revisit the Ising model with three-spin interactions on a triangular 
lattice, called Baxter- Wu (BW) model after its exact solution obtained thirty years ago [9] . 
The motivation arises from the disagreement between results for the exponent z obtained for 
the BW model and for two other models - the two-dimensional four-state Potts model and 
the Ising model with three-spin interactions in one direction - which are known to belong to 
the same universality class. The BW model consists of Ising spins (Si = ±1) on a triangular 
lattice with the Hamiltonian 

H = -J ]T SiSjSk, (1) 

<i,j,k> 

where the sum involves the product of the spins on each triangle of the lattice. As the Ising 
model on a square lattice this model undergoes a phase transition at the critical temperature 
K c = J/kT c = 0.51n(l + \^2). The ground-state is fourfold degenerate (see Fig. 1) and the 
relevant symmetry is semi-global (the Hamiltonian is invariant under reversal of all of the 
spins belonging to two of three sublatices into which the original lattice can be decomposed). 
The BW model exhibits the same static critical exponents (a — v — 2/3 and r\ = 1/4) [10] as 
those of the two-dimensional four-state Potts model and of the Ising model with three-spin 
interactions (IMTSI) in one direction [11]. However, finite-size studies of the BW model have 
good convergence, without running into the known difficulties encountered in the studies of 
these other models. Several authors [12] have attributed that difference to a dangerous 
marginal operator present in the 4-state Potts and IMTSI models which does not have a 
counterpart in the BW model. 

Although the static critical behavior of these systems is well understood, only a few 
results are known regarding dynamical properties. For instance, much work it is needed to 
unravel the interesting question of a possible extended universality, i. e., if these three models 
(BW, four-state Potts and the IMTSI) exhibit the same exponent z. Previous estimates of 
the critical exponent z for the four-state Potts model are scattered between 2 and 4 [13]. The 
most recent results however suggest values close to 2.3 [14,15]. Short-time simulations for 
the BW model were first performed by Santos and Figueiredo [16], who found z = 2.07(1), 
a very small value when compared to the values obtained for the same exponents of the 
four-state Potts and IMTSI models [14,15]. Moreover, they stress that this value of z is very 
close to the dynamic exponent of the the two-dimensional (2D) Ising model when Glauber 
dynamics is used. For the new dynamic exponent 9 related to the critical initial slip of 
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the order parameter, situation is even worse. Okano et al have conjectured [4] that for the 
four-state Potts model the exponent 9 should be negative and close to zero. This result 
was only recently confirmed [18]. For the IMTSI model, two independent estimates were 
obtained [15,19] and fall in the range —0.04 < 9 < 0.07 , whereas estimates for the BW 
model are not available. 

Motivated by our recent results for the IMTSI model [15], we decided to investigate the 
short-time critical dynamics of the B W model in order to check the value of the dynamic 
exponent z and the utility of the exponent 9 in detecting the presence of a marginal operator 
[5,6,20,27]. In the next section we calculate the exponent 9 studying the evolution of the or- 
der parameter for different values of the initial sublattice magnetization m = (3/iV) J2i S'i(O) 
where the index % runs over only one kind of sites of the triangular lattice (circle, square 
or star points in Fig. 1) and iV is the number of spins. In addition, we estimate the same 
exponent calculating the time correlation of the order parameter [21]. In section III we 
summarize our results for the dynamic critical exponent z. In section IV we list estimates 
for the correlation length exponent v and for the index (3 of the magnetization. Finally, in 
section V we present our conclusions. 



Janssen et al. [1] have shown on the basis of the renormalization group theory that, at an 
early stage of the time evolution, quantities such as the magnetization and its moments follow 
a power law when the parameters of a dynamical system are adjusted to their critical values. 
The exponent 9 which characterizes the time evolution of the magnetization is independent 
of the known set of static critical exponents and even from the dynamic exponent z. However, 
like other exponents 9 seems to depend on general properties such as dimensionality of the 
spin and either symmetry of the Hamiltonian or the evolution rules. Such universal behavior 
was found in two and three-dimensional Ising models [4,22], two-dimensional three and four- 
state Potts models [14,23] as well as in irreversible models with synchronous and continuous 
dynamics [21,24,25]. In addition, nonuniversal behavior was detected for the exponent 9 
[5,6,26,27] although always accompanied by a marginal operator. 

In this paper we estimate the dynamic exponent 9 by two different techniques. First we 
quench the system from a high temperature to T c and follow the nonequilibrium relaxation 
of various realizations of the process with the same initial sublattice magnetization m . In 
Fig. 2 we show the behavior of the order parameter when m = 0.01. As proposed by 
Janssen et al [1] it follows the power law 



but here 9 is negative. Figure 3 exhibits the behavior of the exponent 9 for different values 
of mo, as well as a linear fit which leads to the value 9 = —0.186 ± 0.002. 

The second technique consists of studying the time correlation of the magnetization 



which was shown [21] to behave as t 6 when < M(0) >= 0. In Fig. 4 we show the time 
dependence of the time correlation C. We stress that this approach was shown to be valid 



II. THE DYNAMIC EXPONENT 9 



M(t) ~ mot , 



(2) 




(3) 
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when the evolution matrix has up-down symmetry which is not valid in the present case. 
However, it is possible to show [28] that the Z(2) ® Z{2) symmetry of the BW Hamiltonian 
is sufficient to preserve the above mentioned result. The value of 9 obtained with this 
technique (9 = —0.185 ± 0.001) corroborates the value obtained directly from the behavior 
of the magnetization. We stress that the time correlation method [21] is numerically more 
convenient since we do not need to work with specially prepared initial configurations (with 
fixed mo) nor it is necessary to take the delicate limit mo — > 0. 

The results obtained with the two techniques are essentialy the same and fall in the 
range -0.188 < 9 < -0.184. They are completely different from the value 9 = -0.03 ±0.01 
encountered for the IMTSI model [15] and from 9 = -0.047 ± 0.0033 [18] obtained recently 
for the four-state Potts model. Thus, based on our estimates for the exponent 9 of the BW 
model, we conclude that its dynamic behavior at short-time is completely different from 
that of the IMTSI and four-state Potts models albeit both models have the same static 
exponents as the BW model. Thus, at least in this case, the exponent 9 seems to be able 
to capture the subtle difference between models that have the same leading exponents but 
not the same corrections to scaling. 



III. THE DYNAMIC CRITICAL EXPONENT Z 



Santos and Figueiredo [16] obtained the exponent z of the BW model by means of the 
second order cumulant U 2 [8] 

< M 2 > 

starting from completely ordered samples (mo = 1)- As claimed by Zheng, U 2 should behave 
at the critical temperature as a power law 

U 2 oc t d '\ (5) 

which would allow us to estimate z without using any other critical exponents. Santos 
and Figueiredo found z = 2.07(1) for m = 1 (all sub-lattices initially with spins up) and 
1.96(2) for mo = —1/3, i.e., spins are initially up in one sub-lattice and down in the other 
two. These estimates, obtained with Glauber dynamics [17], are very different from results 
recently obtained for the four-state Potts model [14] and for the Ising model with three spin 
interactions in one direction [15] which both yield values close to 2.3. 

This kind of disagreement has already appeared in the literature when the second cumu- 
lant was used to obtain estimates for the exponent z for the 3— state Potts model [8] and 
the IMTSI model [15]. Thus, in order to check if the BW model (which also belongs to the 
4— state Potts model universality class) exhibits indeed a different value for the dynamic 
exponent z, we have calculated this exponent using three other techniques. To begin with, 
we estimated z by collapsing the generalized fourth-order Binder cumulant 

< M 4 > 
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which should only depend on the ratio t/L z when T = T c , according to scaling laws valid 
in the beginning of the evolution. In Fig. 5 we show the Binder cumulant as a function of 
the time for two different lattices L = 96 and L' = 192, both starting from ordered states 
(m = 1). We also plot in Fig. 5 the cumulant for the larger lattice against the re-scaled 
time t' = 2 z t where z is defined via the relation 

U 4 (2 z t,0,L') = U 4 (t,0,L). (7) 

After interpolating the re-scaled curve for the lattice V , the value of z is selected using 
the x 2 test [29] yielding the best estimate 2.285(10). We applied the same analysis to the 
cumulant [8] 

< M 2 > 

U= l2 -1, (8) 
|M| 2 V 1 

using (mo = 1) as shown in Fig. 6. The result in this case is z = 2.29(1). 

Finally, we used mixed initial conditions in order to obtain the function F 2 (t), given by 

which was shown [14] to behave as t d l z . This approach was recently proposed in order to 
obtain the exponent z with high precision. Actually, curves of F 2 are straight lines (see Fig. 
7) which do not depend on the lattice size (one does not need to rescale the time). Our 
estimate for z using Eq.(9) and taking values at each 5 steps is z = 2.294(6). This value is 
in complete agreement with the result z = 2.290(3) recently obtained [14] for the four-state 
Potts model and corroborates the above results obtained by different procedures. 

It is worth stressing that all the results shown until now were obtained with heat-bath 
updating whereas Santos and Figueiredo [16] performed simulations using Glauber dynamics. 
In order to rule out any influence of the dynamics on the results we repeated our simulations 
using Glauber updating and verified that our result remain essentially the same. In addition, 
we have used heat-bath updating to calculate the time dependence of the second order 
cumulant (Eq.(4)) where the average is always obtained from initially ordered samples (mo = 
1). Our best estimate for z is 2.03 ± 0.01, which agrees reasonably well with the result by 
Santos and Figueiredo, indicating that this severe disagreement stems from a difficulty from 
the scaling for the second cumulant method. 

IV. THE THERMAL EXPONENTS v AND (3 

Taking the logarithmic derivative of the magnetization [8] 

M(t, t) = r"/ M M(l, t^r) (10) 
with respect to r = (T — T c )/T c , we obtain 

d T ]nM(t,r) = t x l vz d T , ln(F(r') | T , =0 , (11) 
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which follows a power law that does not depend on L, where F(t') is a scaling function. Thus, 
a log- log plot of that derivative versus time gives 1/vz (see Fig. 8). Once we have obtained 
the dynamical exponent z, the thermal exponent v may be determined. Our estimate for 
the exponent v in the short-time analysis is 0.65 ± 0.01 which should be compared to the 
exact value 2/3 from Ref. [30]. 

Next, we follow the decay of the sublattice magnetization in initially ordered samples 
(m = 1). The scaling law 

M oc t- p/vz (12) 

allows us to obtain fi/vz which in turn leads to the exponent (5 by using the previous result 
for the product vz. We show a log-log plot of M versus t in Fig. 9. A linear fit gives the 
value 0.0527 ± 0.0005 for the ratio (3/vz leading to the estimate 

P = 0.079 ± 0.002. (13) 

We stress that our result for (3 is in fair agreement with the exact result 1/12 [9]. 

Before closing this section we present another estimate for the exponents f3 and v obtained 
from caling of the so-called Baxter order parameter Mb given by 

M B = ]T SiSjS k , (14) 

<ijk> 

where the sum is done over all triangles of the lattice and the indices i, j and k denote sites 
composing an elementary triangle (plaquette). Barber and Baxter [31] have shown that 
Mb is zero for temperatures below T c and different from zero above T c . In this sense Mb 
functions as a disorder parameter [32] which goes to zero when T — > T c from above. However, 
this new quantity has the same critical exponent as the sublattice magnetization. Thus, the 
exponent (3 can be obtained following the evolution of the Baxter's order parameter. The 
advantage is that we do not need to divide the system into sub-lattices. In Fig. 10 we show 
the decay of the Baxter order parameter as well as the value of the ratio (3/vz. In order 
to eliminate the product vz from our result we start again from ordered samples (M B = 1) 
and use Eq.(ll) replacing M by M B to obtain 1/vz. Figure 11 shows the evolution of the 
derivative. Our new estimate for 1/vz is 0.651 ±0.006, yelding for v the value 0.67 ±0.01 in 
very good agreement with the exact result 2/3 [9]. The new estimate for f3 is 0.080 ± 0.002 
which is in fair agreement with our previous result and should be compared to the pertinent 
result 1/12. 

V. CONCLUSIONS 

We revisited the Baxter- Wu model using short-time Monte Carlo simulations. Our best 
estimate for the dynamic critical exponent is z = 2.294 ± 0.006. It was obtained by a 
recently proposed technique [14] that combines simulations performed with different initial 
conditions. This result is compatible with recently obtained estimates for z in the case of the 
two-dimensional four-state Potts model [14] as well as for the two-dimensional Ising model 
with three-spin interactions along one direction [15]. The result for z was confirmed by two 
other techniques based on the collapse of cumulants. In addition, using scaling relations 
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for the order parameter and its derivative with respect to the reduced temperature r, we 
were also able to estimate the static critical exponents (3 and v. We have also obtained the 
exponents (3 and v using the Baxter's order parameter. The results for these exponents are 
in good agreement with the exact values [9]. However, our estimates for the new dynamic 
exponent 9 of the Baxter- Wu model are completely different from the values found for that 
index when studying the two-dimensional Ising model with three-spin interactions in one 
direction and the four-state Potts model. In this sense, the exponent 9 seems to be able to 
capture subtle difference between models which exhibit the same leading critical exponents 
but different corrections to scaling. 
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Figure Captions 



Figure 1. The Baxter- Wu model is defined on a triangular lattice where Ising spins 
(Si = ±1 ) are coupled by nearest neighbor three-spin interactions. The lattice can be 
divided in three sub-lattices, characterized in the figure by circle, square and star points. 
Each symbol represents a class of spins which does not interact with each other. The ground- 
state is fourfold degenerate, corresponding to the three configurations with up spins in two 
sub-lattices and down spins in the other, besides the obvious configuration of spins up in all 
the sublattices. 

Figure 2. Power law behavior of the magnetization for samples with initial magnetization 
m = 0.01. The error bars were calculated over 10 sets of 60000 samples. 

Figure 3. Dynamic exponent 6 as a function of the initial magnetization m for the 
lattice size L = 60. Each point represents an average over 10 sets of 60000 samples. 

Figure 4. Time correlation of the total magnetization for samples with < M(t = 0) >= 0. 
Error bars were calculated over 10 sets of 60000 samples. 

Figure 5. Binder cumulant U±(t, 0, L) for lattice size L = 96 and V = 192 and initial 
magnetization m = 1. The open circles on the continuous line show the cumulant for the 
bigger lattice re-scaled in time with z obtained from Eq. (7). 

Figure 6. U cumulant for L = 96, V = 192 and initial magnetization mo = 1. The open 
circles correspond to the re-scaled cumulant for the larger lattice with z = 2.29(1). 

Figure 7. Time evolution of F2 for L = 102 with mixed initial conditions [Eq.(9)]. The 
error bars, calculated over 10 sets of 30000 samples, are smaller than the size of the points. 

Figure 8. Logarithmic derivative of the magnetization with respect to the reduced tem- 
perature DM(t) as a function of time. Data were obtained for the lattice size L = 102 using 
temperature variation K c r = 0.001. The slope corresponds to the ratio 1/vz, after Eq. (11). 
The error bars, calculated over 10 sets of 30000 samples, are smaller than the size of the 
points. 

Figure 9. Log-log plot of the time evolution of the magnetization for initially ordered 
samples (m = 1). The error bars (calculated over 5 sets of 1000 samples) are smaller than 
the size of the points. 

Figure 10. Time decay of the Baxter's order parameter Mb- Error bars were calculated 
over 5 sets of 10000 samples. 

Figure 11. Logarithmic derivative of the Baxter's order parameter M B . Data were 
obtained for the lattice size L = 204 using temperature variation K c r = 0.001. Five sets of 
10000 initially ordered samples (m = 1) were used in order to obtain the error bars. 
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